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Abstract 

Given the cost and complication of simulating Stirling convertors, do we really need 
multidimensional modeling when one-dimensional capabilities exist? This paper provides 
a comprehensive description of when and why multidimensional simulation is needed. 


I. Introduction 

P OWER conversion with free-piston Stirling engines 1,2 has the potential to deliver high efficiency, low 
mass solutions for longer and more varied space missions . 3,4 The design of free-piston Stirling engines 
has largely been accomplished utilizing one-dimensional Navier-Stokes solvers. One of the objectives of 
computer analysis is to determine an approximate optimum design prior to building any hardware. A 
one-dimensional analysis design code does not provide all of the geometrical details, but can be used to 
determine overall volume and configuration information, including some of the geometrical details. The one- 
dimensional models can be set up quickly and the computations are fast. And perhaps most importantly, 
design optimizations are easily done in one-dimension. 

This paper will demonstrate that one-dimensional analysis has an important place in Stirling engine 
design, but that the design process really needs multidimensional capability as well to: 

• verify the one-dimensional results, 

• properly simulate inherently three-dimensional turbulence-including transition, (includes understand- 
ing the effects of flow vortices in expansion, compression and other spaces on engine performance), 

• provide empirical heat transfer and friction factor coefficients for complex geometries to aid in hardware 
design and to avoid requiring hardware tests to improve one-dimensional analysis assumptions, 

• integrate all the parts into a CAD system and test for structural and relative motion issues in the 
overall design, 

• assist experimentalists in determining flow and heat transfer physics in regions where measurements 
are difficult or impossible, 

• provide a fluid- structure interaction capability, 

• generate linear reduced order models for controller development 
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• model large, high-power Stirling applications or devices in which the one-dimensional flow assumption 
breaks down (even some relatively low power, low At, designs approach a ” pancake” shape that could 
not be well modeled with one-dimensional codes, 

• identify areas of excessive flow losses due to unintended dead zones, recirculation zones, dissipative 
turbulence and other losses such as: 5,6,7 

1. Inefficient heat exchange and pressure loss in the regenerator, heater and cooler, 

2. Gas spring and working space loss due to hysterisis and turbulence, 

3. Appendix gap losses due to pumping and shuttle effects, 

4. Mixing gas losses from nonuniform temperature and flow distributions perpendicular to primary 
engine flow axis, 

5. Conduction losses from the hot to cold regions 

6. Losses due to combined radiation, conduction and convection in void volumes 

7. And in general, inaccurate loss representations due to use of 1-D flow design codes to account 
for flow and heat transfer through area changes (between components) where phenomena such as 
flow separation and jetting from tubes or slots into a regenerator may occur. 

This paper will highlight why multidimensional Stirling analysis is useful and critical for future designs. 

II. Description of the Problem 

The dual opposed configuration shown in Fig. 1 8,9,10 j s b e i n g developed for multimission use (i.e., for 
use in atmospheres and space), including providing electric power for potential missions such as unmanned 
Mars rovers and deep space missions. 11 



Figure 1. Dual Opposed Stirling Convertors Reduce Vibration (Schreiber) 

Only the Stirling engine part of the convertor (Fig. 2) is simulated multi-dimensionally although one 
could anticipate the entire convertor may one day be prototyped digitally before build-up. 


NAS A/TM— 2005-2 13975 


2 






Figure 2. Stirling Engine Part of Convertor Simulated Only 

An example of the desired full Stirling engine simulation is shown in 
Fig. 3. 12 It is anticipated that full 3D simulations will provide a level of 
geometric and flow detail necessary for further design improvements. For 
example, hardware experiments have shown that large performance gains 
can be made by varying manifolds and heat exchanger designs to improve 
flow distributions in the heat exchangers. 8, 13, 14 

The kinds of flow and geometry that occur in a Stirling engine which 
a modeling technique must address are as follows: 7,6 

1. Oscillating flow which changes the effective conduction, flow friction 
coefficients and heat diffusion. 15 

2. Low mach number flow (no shocks), 

3. Compressible flow due to enclosed varying volumes and heat transfer 
effects, 

4. Laminar, Transitional, and Turbulent flow with Reynolds numbers 
from 100 to 10000 (based on various length scales pertinent to the flow region), 

5. Conjugate heat transfer, thermal dispersion and local thermal nonequi- 
librium in the porous media regenerator, 

6. Micron to Millimeter scale geometry 16 

7. Small separations between sliding solid-fluid interfaces in the seals 
and appendix gap 

8. Deforming Flow Regions Due to Compression and Expansion (i.e., 
due to moving piston boundaries) 

III. Current Stirling Cycle Analysis Approaches 

The techniques of analysis for Stirling engines can be categorized with Martini’s 17, 18 nomenclature. We 
will use his definition of Third Order and discuss the recent multidimensional category. 


Figure 3. Example Transient 
Simulation of Stirling Convertor- 
Colored by Temperature 
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A. Third Order Analysis 

Third order analysis uses control volumes or nodes to directly solve the one-dimensional governing equations. 
Some of the first analysis at this level of fidelity was by Finkelstein, 19 Urieli, 20 and Berchowitz. 21 

Some other more recent third order analysis codes are: 

• The codes by David Gedeon referred to as GLIMPS 22,23 and Sage 24,25,26 are one-dimensional and solve 
the governing equations implicitly in space and time. The grid includes all time because a periodic 
solution is assumed/forced. Therefore, it is not possible to model transient startup behavior. 

• The linearized harmonic analysis code referred to as HFAST 27 solved a steady-state periodic problem 
in the frequency domain. Again, transient behavior is not modeled. 

• The code by Martini engineering 28 was never validated but claimed transient modeling capability. It 
is not clear that rigorous governing equations are being solved since it appears many simplifications 
based on experimental correlations are used. 

• Another unvalidated but interesting code by Renfroe 29 attempted a one-dimensional analysis using 
explicit Runge-Kutta time stepping and a Newton solver to solve the nonlinear equations. 

• Some other recent one-dimensional solvers are LASER, 30 DeltaE, 31 ARCOPTR, 32 and REGEN3.1. 33 

• Finally, the Stirling Dynamic Model (SDM), 34 uses a one-dimensional analogy of an entire Stirling 
convertor by linking together representative elements within the Simplorer(TM) commercial software 
package by Ansoft Corp. This tool enables approximate whole convertor dynamic analysis. Recent 
work attempts to incorporate thermodynamics via David Gedeon’s Sage code described earlier. 

These approaches will continue to play an important role in Stirling analysis, even as the multidimensional 
analysis becomes practical to use. A full analysis package incorporating both third and fourth order analysis 
(defined below) will be presented in Section VI. 

B. Fourth Order (Multi-Dimensional) Analysis 

At this level of analysis, relatively little has been completed because the third order analysis is faster and for 
the most part has been an adequate engineering tool. However, to improve efficiency further (to understand 
and reduce losses) it will likely require a better understanding of the actual flow physics and heat transfer 
throughout the engine. It is impractical to measure many of these features in practical Stirling devices. 
Some of these features can be investigated in large scale test modules designed to simulate certain Stirling 
like processes (such as those being investigated at the University of Minnesota. However, CFD seems best 
suited for investigation of the details of the fluid and heat transfer physics in a real Stirling device. Some 
multi-physics analysis tools are listed below. 

1. Modified Computer Aided Simulation of Turbulence (CAST) 

The modified CAST code 35 is based upon the Semi-Implicit Method for Pressure-Linked Equations SIM- 
PLE 36 method but is restricted to two dimensions. 37 It was modified to include oscillatory boundary condi- 
tions and conjugate heat transfer. It has been used to model Stirling components but has not been extended 
to a whole engine simulation tool. 

A pressure-splitting technique was added 38 to reduce the computational requirements. It was based on 
separating the thermodynamic and hydrodynamic pressures so that these widely varying scales could be 
solved with less round-off error and better efficiency. 
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2. CFD-ACE 


This commercial code has been used to model a two-dimensional representative Stirling engine. 39,40 It is 
also based upon the SIMPLE technique. The regenerator is not currently modeled correctly since thermal 
equilibrium is assumed between the gas and solid. 

This finite volume code can utilize both structured and unstructured grids. 

3. Fluent 

This commercial code is also based upon the SIMPLE method (and PISO method for high speed flows). It 
currently has similar regenerator modeling limitations in that it is designed for non-oscillating flows. It does 
however have a sliding interface that could be used for appendix gap modeling on parallel computers. It 
is being used by several commercial manufacturers for this purpose (references are proprietary). It is also 
finite volume based and can utilize both structured and unstructured grids. 

I STAR-CD 

The Simulation of Turbulent Flow in Arbitrary Regions (STAR) 41 code also uses SIMPLE and PISO methods. 
Its companion product (STAR-HPC) is the parallel computer version. It also has sliding interfaces and 
deforming mesh capability. It has been used in the related field of internal combustion engine piston modeling, 
and some Stirling engines have been modeled with it. 42 

5. CFX 

This code also uses SIMPLE and PISO methods on unstructured grids. It also has sliding interfaces imple- 
mented, but no Stirling engine modeling with this software has been publicly published. 43 

6. Others 

While there are other in-house codes, they are usually limited to modeling only specific regions of the Stirling 
engine such as the regenerator, or the displacer. 

IV. Recent Whole Engine Modeling 

The use of two-dimensional CFD models can significantly extend the capabilities, compared to third-order 
analysis, for the more detailed analysis of the complex heat transfer and gas dynamical processes which occur 
in the internal gas circuit. 44 More recently, full 3-D calculations have been performed with a commercial 
code. 45 The temperature results are similar to the second order method results. 46 The multidimensionally 
computed power, however, was about half of the second order prediction. Moreover, along the axis of the 
compression space it was found the change of the temperature of the working gas was quite different from 
harmonic in time. 

Zhang 47 claims success with modeling a 3D free-piston pseudo- Stirling engine over a 3 month run-time. 
More recently, a two-dimensional axisymmetric simulation of a full engine has been demonstrated and vali- 
dated within the observed experimental results taken from two engines. 12 Moreover, it is possible to simulate 
an entire engine cycle in less than one hour by utilizing modern parallel computer architectures. 

In short, the capability does now exist to perform whole engine simulations in considerably less time than 
previously expected. 
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A. Regenerator Modeling 



A very important specific area of modeling difficulty is the re- 
generator. As shown in Fig. 4, since the regenerator (depend- 
ing upon one’s definition of effective) has roughly 3 to 40 times 
more effective heat transfer than the heater, 48 any inefficiency 
of the regenerator represents a significant loss for the entire 
Stirling engine. Hence, any numerical losses or inaccuracies in 
this region will disproportionately influence the entire Stirling 
simulation. 

1 . Manifest 

The multidimensional code referred to as Manifest 16 solved the 
porous medium model equations shown below: 


Time 


Figure 4. Importance of Properly Modeling 
the Regenerator for System Studies (Urieli) 


dM , J d ( G k \ _ n 
dr ^ (3 dy k \ J ) ~ U 

t §*) + (? §») + «.*»*) a* 

W + iA (J/ 3 ”* + W - )(r. - 13V)" + -Q = 0 
£(\t s ) + q = 0 


(1) 

where M = p, G = /3pV , E = pe, using the Beam and Warm- 
ing 49 implicit time stepping approach. 

This code was written in curvilinear coordinates and seemed 
to perform well for low Reynolds number jets impinging on the 
matrix. It was not developed further because of excessive computational time. However, recent progress in 
improvements to efficiency suggest revisiting this approach. 12,50 




(c) Photographed Geometry 


Figure 5. Regenerator Geometry 


Another important issue is the geometrical shape of the matrix in the regenerator. As shown in Fig. 5, 
most regenerator models don’t assume a precise geometrical shape for the elements of the regenerator. 
However, as reported in Park, 51 a simple-to-fabricate woven mesh, consisting of bonded laminates of two- 
dimensional plain-weave conductive screens can be manufactured to have a wide range of porosity and a 
highly anisotropic thermal conductivity vector. In addition to providing superior performance in many 
cases, 52 the regular geometry greatly simplifies the analysis. 
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Clearly, the shape of the regenerator has an important impact on the overall system design and attempts 
have been made in the past to multidimensionally simulate the different mediums. 

V. Need for Multi-Dimensional Analysis Outlined 

Recent accomplishments in multi-dimensional simulation mentioned above suggest fast solutions are 
possible and the value of these results lie in the interplay between supporting current one-dimensional 
results and providing information previously unattainable. For example, the following areas are ripe for 
multidimensional analysis. 

A. Turbulence Modeling 

Since turbulence is random and inherently three-dimensional, a quasi-steady harmonic solution is capable of 
fully modeling its effect on pressure drop and heat transfer. A full three-dimensional Large Eddy Simulation 
(LES) simulation would be capable and could be utilized to determine some of the heat transfer coefficients 
and friction factors prior to any hardware experiments. Modeling transition effects on these coefficients 
appears possible by using recently developed numerical techniques. 50 This information could then be utilized 
by the one-dimensional solver to simulate engine components that have not been built before. 

B. Verify One-Dimensional Results 

The less expensive one-dimensional results are fairly accurate in most cases, so long as the empirical coeffi- 
cients used are accurate. After a preliminary design is optimized with a one-dimensional code, a check on 
the final design with a full three-dimensional simulation would eliminate the expense and risk of building an 
engine that did not perform as expected when using only one-dimensional assumptions. 

C. Numerical Empirical Coefficient Determination 

Similar to the advantages of proper turbulence modeling, certain empirical heat transfer and flow friction 
factors could be obtained with multidimensional analysis when the geometry is sufficiently deviant from one- 
dimensional flow physics. For example, the flat head heater design by Infinia (formerly STC) was initially 
modeled using a one-dimensional solver but due to the complex geometry (Fig. 6) the error was over 20% 
until the empirical terms were adjusted. 53 Other examples include high power applications with designs that 
are radially wide compared to the axial length. 

D. Part Integration 

A conceptual one-dimensional design still does not provide adequate information for drawing actual hardware. 
An important aspect of a good design is not only the engine efficiency, but also it must be light-weight, 
reliable, and low cost. The ability to test several orientations of the parts before building them is another 
advantage of fourth order modeling. For example, in Fig. 7, a preliminary CAD drawing (assembly of parts) 
of a piston configuration and related part integration is shown. 

E. Assisting Experiments 

Experimental uncertainty can be minimized by multi-dimensionally simulating the experiment. Regions 
inside the engine are difficult to measure due to the small area or because the sensors would interfere with 
the flow. For example, knowledge of the temperature range experienced by the flexures throughout a cycle 
is important for determining the reliability of the convertor over long periods of operation in space. A full 
dynamic measurement of temperature distribution throughout a Stirling device is not practical due to sensor 
limitations and the reasons just mentioned. However, it is a fairly simple matter with multidimensional 
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S. Qui, STC, IECEC 2004 


Figure 6. Flat Head Heater Not One-Dimensional 


simulations to access this data. In Fig. 8, a thermocouple cannot easily be located on the cold wall of the 
engine and so precise temperature information for the operating conditions of the engine are not known. 
However, a three-dimensional simulation has been used to identify correlations for the coldwall temperature 
inside a cooling jacket. 

The temperature and velocity distribution throughout the engine is important for determining the kinds 
of materials and joining techniques. This information is also most easily determined numerically. An example 
of this for a simple linear alternator piston problem at one moment in time is shown in Fig. 9. 

F. Fluid-Structure Interaction 

The fluid forces will bend the structural materials and the structure will in turn transfer heat with the fluid 
and push the fluid. For example, inside a typical engine, the radiation shields in the displacer are very thin 
and any flexing will change the convection path of the neighboring fluid. Similarly, the flexures will mix the 
Helium gas changing the effective temperature on the flexure. These, and other interactions are inherently 
multi-dimensional. 

G. Dimensionality of Losses 

The primary loss mechanisms may be grouped as mixing, viscosity, and irreversible heat transfer. Mixing is 
inherently multidimensional. Losses due to viscosity are often found from friction factors determined from 
experiment. And the losses from heat transfer depend upon heat transfer (film) coefficients that are also 
experimentally derived. A multidimensional analysis does not require the experimental data. 
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Figure 7. Fitting Parts Together 


H. Generating Linear Models for Controls 

Computational Fluid Dynamics based linear modeling provides a small perturbation model that can be used 
for control applications and real-time simulations. 54 By producing a steady-harmonic solution with a multi- 
dimensional analysis, small changes in the operating conditions of the engine will produce small changes in 
performance. This information can be used in developing controllers for dynamic response. 

VI. Design And Integration Analysis Options 

The capability for two-dimensional axisymmetric multidimensional Stirling simulation has recently been 
demonstrated and the solutions take about a week. However, through a careful integration of multiple 
solution strategies, a more complete design analysis may be achieved in less time. For example, in Fig. 10, 
is a road map for achieving a more integrated approach to Stirling engine design. 

A. Multi-D Sage 

Currently Sage 24 is a one-dimensional Navier-Stokes steady- harmonic solver commonly used in industry for 
Stirling engine performance analysis and optimization. It’s main weakness is its one-dimensional formulation. 
The solution approach could be extended to two and three-spatial dimensions for very fast solutions of 
complex engines of all power levels. This would also provide correct boundary conditions for component 
modeling 55 since current one-dimensional analysis does not provide multidimensional velocity direction. 

B. Commercial Transient Code Upgrades 

Current commercial codes utilize old techniques that were designed for steady-state calculations. However, 
the commercial codes are well developed, user friendly, and useful for simulating complicated moving ge- 
ometry such as occurs in Stirling engines. It is important to utilize this resource for near term results and 
modest improvements to their regenerator modeling equations would make them useful for design this year. 

They can provide improved heat transfer and friction factor correlations for one-dimensional Sage cal- 
culations. In addition, gas spring and flexure spring coefficients could be determined multidimensionally 
for improved System Dynamic Modeling. And finally, an as yet untapped ability is the recently developed 
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Figure 8. Numerically Derived Correlations 





(a) Temperature Distribution 


(b) Velocity Distribution 



Figure 9. Numerically Derived Temperature and Velocity Distributions 


Fluid-Structure interaction capability of the commercial codes. This is still unproven and time-consuming, 
but offers the possibility of reducing displacer rubbing and improved material selection. 

C. High-Order Transient Code Development 

A new high-order transient code would be able to efficiently resolve entropy features and turbulence transi- 
tion. A one-dimensional high order transient code could be incorporated into the SDM to provide thermo- 
dynamic as well as one-dimensional system mechanical modeling. And a multidimensional transient code 
would provide improved correlations for Sage. 

D. SDM-CFD Integration 

Stirling Dynamic Modeling was mentioned earlier and recent attempts at merging it with Sage have proved 
successful. However, higher fidelity simulation needs will naturally lead to incorporating multidimensional 
analysis. An integrated CFD-SDM tool could provide the piston dynamic information for the other analysis 
tools. By iterating between these tools we would converge to a fully validated understanding of the Stirling 
design and this is the subject of the next section. 
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Figure 10. Overall Analysis Design 


E. Comprehensive Analysis 

The need for multidimensional analysis is great, but one-dimensional analysis is important for quickly pro- 
viding conceptual designs and optimizations. By combining all these forms of analysis, we could be assured 
that the first hardware prototype will perform as intended. 

As shown in Fig. 11, a logical design path starting with the desired performance and ending with hardware 
may be followed to generate reliable system solutions. First, after specifying a desired mission power and 
operating range, a preliminary Sage (or other one-dimensional solver) analysis could be performed. The 
motion of the piston could then be improved with the SDM analysis. Using those two results, an initial 
set of CAD drawings of a conceptual engine could be quickly set up and tested in either a two-dimensional 
axisymmetric or full three-dimensional mode. If the expected results are observed then it can be assumed 
the empirical coefficients used in Sage and SDM are adequate. Next, optimize the engine in Sage and again 
repeat the SDM analysis. Repeat the 2D and/or 3D analysis to confirm the results. If at any time the one- 
dimensional results do not match the multidimensional results then change the emperical coefficients based 
on the multidimensional results. This process of iterating between one-dimensional and multidimensional 
results insures the hardware will function as intended. 
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Figure 11. Analysis Integration 


VII. Conclusion 

These observations lead to the conclusion that one, two, and three-dimensional modeling should all be 
employed and all three paradigms provide important capabilities that when combined provide a potent com- 
bination of initial design, empirical coefficient adjustment, optimization, and final prototype demonstration 
before the first part is cut. 
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